knitr::opts_chunk$set(echo = TRUE)

cmdstanr::set_cmdstan_path(path = "C:/Users/kueng/.cmdstan/cmdstan-2.35.0")


library(tidyverse)
library(R.utils)
library(wbCorr)
library(readxl)
library(kableExtra)
library(brms)
library(bayesplot)
library(beepr)
library(DHARMa)



source(file.path('Functions', 'ReportModels.R'))
source(file.path('Functions', 'PrettyTables.R'))
source(file.path('Functions', 'ReportMeasures.R'))
source(file.path('Functions', 'PrepareData.R'))
system("shutdown /a")
## [1] 1116
# Set options for analysis
use_mi = FALSE
shutdown = TRUE
report_ordinal = FALSE

iterations = 2000
#iterations = 10000

options(
  dplyr.print_max = 100, 
  brms.backend = 'cmdstan',
  brms.file_refit = ifelse(use_mi, 'never', 'on_change'),
  brms.file_refit = 'on_change',
  #brms.file_refit = 'always',
  error = function() beepr::beep(sound = 5)
)
df <- openxlsx::read.xlsx(file.path('long.xlsx'))
df_original <- df

# Importantly, we do not recode pushing here. 
df_double <- prepare_data(df, recode_pushing = FALSE, use_mi = use_mi)[[1]]

Constructing scales reshaping data (4field) centering data within and between

Modelling

# For indistinguishable Dyads
model_rows_fixed <- c(
    'Intercept', 
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb'
  )


model_rows_fixed_ordinal <- c(
  model_rows_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rows_fixed[2:length(model_rows_fixed)]
)

model_rows_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)

model_rows_random_ordinal <- c(model_rows_random,'disc')
# For indistinguishable Dyads
model_rownames_fixed <- c(
    'Intercept', 
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'Daily received persuasion target -> target', 
    'Daily received persuasion target -> agent', 
    'Daily received pressure target -> target', 
    'Daily received pressure target -> agent', 
    'Daily received pushing target -> target', 
    'Daily received pushing target -> agent', 
    'Day', 
    'Daily weartime',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'Mean received persuasion target -> target',
    'Mean received persuasion target -> agent',
    'Mean received pressure target -> target',
    'Mean received pressure target -> agent',
    'Mean received pushing target -> target',
    'Mean received pushing target -> agent',
    'Mean weartime'
  )


model_rownames_fixed_ordinal <- c(
  model_rownames_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed[2:length(model_rownames_fixed)]
)

model_rownames_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(Daily received persuasion target -> target)', 
  'sd(Daily received persuasion target -> agent)', 
  'sd(Daily received pressure target -> target)', 
  'sd(Daily received pressure target -> agent)', 
  'sd(Daily received pushing target -> target)', 
  'sd(Daily received pushing target -> agent)', 
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]',
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)

model_rownames_random_ordinal <- c(model_rownames_random,'disc')
# For indistinguishable Dyads
model_rownames_fixed <- c(
    "Intercept", 
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime"
  )


model_rownames_fixed_ordinal <- c(
  model_rownames_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed[2:length(model_rownames_fixed)]
)

model_rownames_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)

model_rownames_random_ordinal <- c(model_rownames_random,'disc')
rows_to_pack <- list(
  "Within-Person Effects" = c(2,9),
  "Between-Person Effects" = c(10,16),
  "Random Effects" = c(17, 23), 
  "Additional Parameters" = c(24,30)
  )


rows_to_pack_ordinal <- list(
  "Intercepts" = c(1,6),
  "Within-Person Effects" = c(2+5,9+5),
  "Between-Person Effects" = c(10+5,16+5),
  "Random Effects" = c(17+5, 23+5), 
  "Additional Parameters" = c(24+5,30+6)
  )
HURDLE MODELS
# For indistinguishable Dyads
model_rows_fixed_hu <- c(
    'Intercept', 
    'hu_Intercept',
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb',
    
    # HURDLE MODEL
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'hu_persuasion_self_cw', 
    'hu_persuasion_partner_cw', 
    'hu_pressure_self_cw', 
    'hu_pressure_partner_cw', 
    'hu_pushing_self_cw', 
    'hu_pushing_partner_cw', 
    'hu_day', 
    'hu_weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'hu_persuasion_self_cb',
    'hu_persuasion_partner_cb',
    'hu_pressure_self_cb',
    'hu_pressure_partner_cb',
    'hu_pushing_self_cb',
    'hu_pushing_partner_cb',
    'hu_weartime_self_cb'
  )

model_rows_random_hu <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(hu_Intercept)',
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # HURDLE
  'sd(hu_persuasion_self_cw)',
  'sd(hu_persuasion_partner_cw)',
  'sd(hu_pressure_self_cw)',
  'sd(hu_pressure_partner_cw)',
  'sd(hu_pushing_self_cw)',
  'sd(hu_pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)
# For indistinguishable Dyads
model_rownames_fixed_hu <- c(
    "Intercept", 
    "Hurdle Intercept",
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime",
    
    # HURDLE
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Hu Daily persuasion experienced", 
    "Hu Daily persuasion utilized (partner's view)", # OR partner received
    "Hu Daily pressure experienced", 
    "Hu Daily pressure utilized (partner's view)", 
    "Hu Daily pushing experienced", 
    "Hu Daily pushing utilized (partner's view)", 
    "Hu Day", 
    "Hu Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Hu Mean persuasion experienced", 
    "Hu Mean persuasion utilized (partner's view)", 
    "Hu Mean pressure experienced", 
    "Hu Mean pressure utilized (partner's view)", 
    "Hu Mean pushing experienced", 
    "Hu Mean pushing utilized (partner's view)", 
    "Hu Mean weartime"
  )


model_rownames_random_hu <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(Hurdle Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  
  # Hurdle
  "sd(Hu Daily persuasion experienced)", 
  "sd(Hu Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Hu Daily pressure experienced)", 
  "sd(Hu Daily pressure utilized (partner's view))", 
  "sd(Hu Daily pushing experienced)", 
  "sd(Hu Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)
rows_to_pack_hu <- list(
  "Conditional Within-Person Effects" = c(3,10),
  "Conditional Between-Person Effects" = c(11,17),
  
  "Hurdle Within-Person Effects" = c(18,25),
  "Hurdle Between-Person Effects" = c(26,32),
  
  "Random Effects" = c(33, 46), 
  "Additional Parameters" = c(47,53)
  )

Subjective MVPA

range(df_double$pa_sub, na.rm = T) 
## [1]   0 720
hist(df_double$pa_sub, breaks = 40) 

hist(log(df_double$pa_sub+00000000001), breaks = 40)

formula <- bf(
  pa_sub ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID),
  
  hu = ~ persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID) 
  #, autocor = ~ cosy(gr = coupleID)
)

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b"),
  brms::set_prior("normal(4, 2)", class = "Intercept"), # for non-zero PA
  brms::set_prior("normal(6, 2.5)", class = "Intercept", dpar = 'hu'), # hurdle part
  
  brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0), # Random effects SD for couple
  
  #brms::set_prior("normal(0, 0.3)", class = "ar", lb = -1, ub = 1), # Autocorrelation prior
  #brms::set_prior("normal(0, 0.3)", class = "ma", lb = -1, ub = 1), # moving average prior
  #brms::set_prior("normal(0, 0.3)", class = "cosy", lb = -1, ub = 1), # compound symmetry prior
  
  #brms::set_prior("normal(0.1,2)", class = "sderr", lb = 0),
  brms::set_prior("normal(0.1,2)", class = "sigma", lb = 0) # Residual SD for lognormal
)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", "pa_sub_hu")
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(pa_sub, log_pp_check = TRUE)
## 
## Divergences:
## 1 of 4000 iterations ended with a divergence (0.025%).
## Try increasing 'adapt_delta' to remove the divergences.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 8 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 4000 by 1342 log-likelihood matrix.
## 
##          Estimate    SE
## elpd_loo  -5682.3  78.5
## p_loo       149.7   6.6
## looic     11364.5 157.0
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.4, 1.5]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     1334  99.4%   212     
##    (0.7, 1]   (bad)         8   0.6%   <NA>    
##    (1, Inf)   (very bad)    0   0.0%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(pa_sub, integer = TRUE)

## Object of Class DHARMa with simulated residuals based on 1000 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help. 
##  
## Scaled residual values: 0.71 0.6028384 0.555 0.1307198 0.2881821 0.894 0.797 0.395 0.958 0.0306795 0.698 0.1053606 0.727 0.376 0.681 0.03561442 0.995 0.995 0.657 0.722 ...
summarize_brms(
  pa_sub, 
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T,
  include_p_direction = T) %>%
  print_df(rows_to_pack = rows_to_pack_hu
)
## Warning: There were 1 divergent transitions after warmup. Increasing
## adapt_delta above 0.8 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Warning in summarize_brms(pa_sub, model_rows_fixed = model_rows_fixed_hu, :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_direction
Intercept 55.11* 45.63 65.99 1.000 1576.42 1787.85 1
Hurdle Intercept 0.34* 0.19 0.62 1.004 1252.89 1442.29 1
Conditional Within-Person Effects
Daily persuasion experienced 1.03 0.98 1.09 1.000 2928.89 3037.14 0.851
Daily persuasion utilized (partner’s view) 1.02 0.97 1.08 1.001 3576.67 2915.87 0.767
Daily pressure experienced 0.89 0.77 1.01 1.000 3326.06 2810.17 0.963
Daily pressure utilized (partner’s view) 0.92 0.81 1.03 1.000 4067.42 2553.90 0.905
Daily pushing experienced 0.98 0.92 1.04 1.001 4599.11 2958.78 0.783
Daily pushing utilized (partner’s view) 0.97 0.91 1.04 1.002 4016.85 3191.01 0.832
Day 0.92 0.79 1.07 1.000 5255.04 3004.93 0.856
Daily weartime NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.08 0.73 1.58 1.001 1293.68 1741.41 0.648
Mean persuasion utilized (partner’s view) 1.04 0.70 1.51 1.001 1355.15 1526.34 0.58
Mean pressure experienced 0.87 0.39 1.95 1.003 1705.28 2397.70 0.641
Mean pressure utilized (partner’s view) 0.71 0.33 1.58 1.002 1672.17 2431.12 0.808
Mean pushing experienced 1.21 0.86 1.70 1.004 1472.49 2325.20 0.863
Mean pushing utilized (partner’s view) 1.18 0.83 1.65 1.006 1371.93 1667.13 0.821
Mean weartime NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 0.65* 0.50 0.81 1.000 3546.54 2574.02 1
Hu Daily persuasion utilized (partner’s view) 0.75* 0.58 0.92 1.001 3522.13 2928.12 0.996
Hu Daily pressure experienced 1.57 0.81 3.00 1.001 3024.74 2505.93 0.921
Hu Daily pressure utilized (partner’s view) 0.69 0.09 1.98 1.003 1513.55 1055.27 0.652
Hu Daily pushing experienced 0.67* 0.44 0.95 1.003 3287.86 2136.06 0.988
Hu Daily pushing utilized (partner’s view) 0.58* 0.36 0.84 1.001 2525.82 2527.75 0.998
Hu Day 1.28 0.74 2.20 1.000 5568.49 2819.54 0.808
Hu Daily weartime NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 0.68 0.14 3.48 1.003 1104.99 1687.61 0.699
Hu Mean persuasion utilized (partner’s view) 0.69 0.15 3.45 1.002 1041.53 1569.39 0.686
Hu Mean pressure experienced 4.05 0.58 37.42 1.001 1433.00 2061.96 0.916
Hu Mean pressure utilized (partner’s view) 3.31 0.49 32.67 1.000 1478.13 2275.30 0.882
Hu Mean pushing experienced 1.56 0.43 5.53 1.002 1279.29 2161.38 0.759
Hu Mean pushing utilized (partner’s view) 1.13 0.32 4.21 1.002 1357.02 2107.68 0.574
Hu Mean weartime NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.37 0.28 0.51 1.00 1291.80 1842.12 NA
sd(Hurdle Intercept) 1.47 1.07 2.01 1.00 1160.68 1834.71 NA
sd(Daily persuasion experienced) 0.09 0.04 0.15 1.00 1796.74 1505.69 NA
sd(Daily persuasion utilized (partner’s view)) 0.08 0.01 0.14 1.00 1327.59 997.32 NA
sd(Daily pressure experienced) 0.09 0.00 0.28 1.00 1864.94 1997.76 NA
sd(Daily pressure utilized (partner’s view)) 0.09 0.00 0.24 1.00 1801.07 1812.28 NA
sd(Daily pushing experienced) 0.05 0.00 0.13 1.00 1622.39 2384.41 NA
sd(Daily pushing utilized (partner’s view)) 0.07 0.00 0.16 1.01 1215.91 1036.82 NA
sd(Hu Daily persuasion experienced) 0.26 0.02 0.60 1.01 1046.40 1633.48 NA
sd(Hu Daily persuasion utilized (partner’s view)) 0.25 0.02 0.56 1.00 1412.08 1778.57 NA
sd(Hu Daily pressure experienced) 0.53 0.01 1.76 1.00 1643.73 2015.28 NA
sd(Hu Daily pressure utilized (partner’s view)) 1.33 0.08 3.87 1.01 1074.00 1605.73 NA
sd(Hu Daily pushing experienced) 0.51 0.06 1.07 1.00 886.92 1237.75 NA
sd(Hu Daily pushing utilized (partner’s view)) 0.55 0.04 1.22 1.00 1073.63 1275.68 NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA
sigma 0.66 0.63 0.69 1.00 4335.61 2618.07 NA

Device Based MVPA

Modelling everything at once

range(df_double$pa_obj, na.rm = T) 
## [1]   5.75 971.25
hist(df_double$pa_obj, breaks = 50)

df_double$pa_obj_log <- log(df_double$pa_obj)

hist(df_double$pa_obj_log, breaks = 50)

formula <- bf(
  pa_obj ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + weartime_self_cw + weartime_self_cb +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = ~ cosy(gr = coupleID)
)


prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b"),
  brms::set_prior("normal(4, 2)", class = "Intercept"), # for non-zero PA
  #brms::set_prior("normal(6, 2.5)", class = "Intercept", dpar = 'hu'), # hurdle part
  
  brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0), # Random effects SD for couple
  
  #brms::set_prior("normal(0, 0.3)", class = "ar", lb = -1, ub = 1), # Autocorrelation prior
  #brms::set_prior("normal(0, 0.3)", class = "ma", lb = -1, ub = 1), # moving average prior
  #brms::set_prior("normal(0, 0.3)", class = "cosy", lb = -1, ub = 1), # compound symmetry prior
  
  #brms::set_prior("normal(0.1,2)", class = "sderr", lb = 0),
  brms::set_prior("normal(0.1,2)", class = "sigma", lb = 0) # Residual SD for lognormal
)



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", "pa_obj_lognormal")
)
## Warning: Rows containing NAs were excluded from the model.
pa_obj_log$prior
check_brms(pa_obj_log, log_pp_check = TRUE)
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 2 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 4000 by 1214 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo  -6911.9 39.4
## p_loo        82.6  5.6
## looic     13823.8 78.9
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.4, 2.3]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     1212  99.8%   326     
##    (0.7, 1]   (bad)         2   0.2%   <NA>    
##    (1, Inf)   (very bad)    0   0.0%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(pa_obj_log, integer = TRUE)
## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## Object of Class DHARMa with simulated residuals based on 1000 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help. 
##  
## Scaled residual values: 0.465 0.616 0.335 0.542 0.67 0.891 0.147 0.898 0.886 0.227 0.71 0.157 0.281 0.286 0.585 0.03 0.815 0.914 0.621 0.617 ...
summarize_brms(
  pa_obj_log, 
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T,
  include_p_direction = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Warning in summarize_brms(pa_obj_log, model_rows_fixed = model_rows_fixed, :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_direction
Intercept 121.95* 107.68 137.86 1.003 954.02 1803.56 1
Within-Person Effects
Daily persuasion experienced 1.01 0.98 1.06 1.002 3255.58 2917.73 0.777
Daily persuasion utilized (partner’s view) 1.02 0.98 1.06 1.000 3409.40 3307.80 0.864
Daily pressure experienced 0.96 0.86 1.06 1.000 3972.57 2940.67 0.811
Daily pressure utilized (partner’s view) 0.94 0.86 1.02 1.002 4676.43 3139.85 0.924
Daily pushing experienced 1.02 0.96 1.09 1.000 4166.74 3316.18 0.79
Daily pushing utilized (partner’s view) 1.01 0.96 1.07 1.003 3647.98 2505.52 0.722
Day 0.99 0.88 1.11 1.001 6757.85 3105.93 0.57
Daily weartime 1.00* 1.00 1.00 1.000 4032.91 2641.61 0.998
Between-Person Effects
Mean persuasion experienced 1.08 0.78 1.50 1.002 923.66 1805.89 0.702
Mean persuasion utilized (partner’s view) 0.89 0.64 1.24 1.003 890.66 1728.13 0.763
Mean pressure experienced 0.77 0.55 1.11 1.001 1478.18 2490.34 0.926
Mean pressure utilized (partner’s view) 0.99 0.71 1.39 1.001 1384.91 2505.65 0.511
Mean pushing experienced 1.09 0.84 1.41 1.002 1018.94 1859.42 0.751
Mean pushing utilized (partner’s view) 1.17 0.90 1.52 1.002 991.68 1746.67 0.888
Mean weartime 1.00 1.00 1.00 1.000 4600.73 3223.11 0.914
Random Effects
sd(Intercept) 0.32 0.25 0.43 1.00 1181.51 1929.75 NA
sd(Daily persuasion experienced) 0.06 0.02 0.11 1.00 1610.76 1611.22 NA
sd(Daily persuasion utilized (partner’s view)) 0.05 0.01 0.11 1.00 1048.24 1121.56 NA
sd(Daily pressure experienced) 0.09 0.00 0.24 1.00 1209.89 1384.88 NA
sd(Daily pressure utilized (partner’s view)) 0.05 0.00 0.15 1.00 2675.11 2349.13 NA
sd(Daily pushing experienced) 0.09 0.01 0.17 1.00 915.18 870.97 NA
sd(Daily pushing utilized (partner’s view)) 0.04 0.00 0.11 1.00 1257.92 1240.76 NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA
sigma 0.54 0.52 0.56 1.00 5268.84 2796.30 NA

Affect

range(df_double$aff, na.rm = T) 
## [1] 1 6
hist(df_double$aff, breaks = 15)

formula <- bf(
  aff ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = ~ cosy(gr = coupleID)
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b"),
  brms::set_prior("normal(3, 2)", class = "Intercept", lb=1, ub=6), # range of the outcome scale
  
  brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0), # Random effects SD for couple
  
  #brms::set_prior("normal(0, 0.3)", class = "ar", lb = -1, ub = 1), # Autocorrelation prior
  #brms::set_prior("normal(0, 0.3)", class = "ma", lb = -1, ub = 1), # moving average prior
  #brms::set_prior("normal(0, 0.3)", class = "cosy", lb = -1, ub = 1), # compound symmetry prior
  
  #brms::set_prior("normal(0.1,2)", class = "sderr", lb = 0),
  brms::set_prior("normal(0.1,2)", class = "sigma", lb = 0) 
)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::skew_normal(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", "mood_skew")
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(mood, log_pp_check = FALSE)
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 50 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 4000 by 1342 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo  -1517.7 35.3
## p_loo        42.9  7.1
## looic      3035.4 70.5
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.5, 1.8]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     1292  96.3%   1129    
##    (0.7, 1]   (bad)        41   3.1%   <NA>    
##    (1, Inf)   (very bad)    9   0.7%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(mood, integer = FALSE)
## qu = 0.5, log(sigma) = -2.401691 : outer Newton did not converge fully.

## Object of Class DHARMa with simulated residuals based on 1000 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help. 
##  
## Scaled residual values: 0.397 0.357 0.659 0.429 0.924 0.691 0.918 0.658 0.683 0.006 0.936 0.815 0.263 0.008 0.442 0.001 0.483 0.885 0.881 0.873 ...
summarize_brms(
  mood, 
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = F,
  include_p_direction = T) %>%
  print_df(rows_to_pack = rows_to_pack)
b l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_direction
Intercept 4.90* 4.80 5.01 1.006 813.62 1375.64 1
Within-Person Effects
Daily persuasion experienced -0.01 -0.04 0.01 1.002 4070.38 3243.64 0.803
Daily persuasion utilized (partner’s view) 0.00 -0.03 0.02 1.000 4386.28 3089.63 0.635
Daily pressure experienced -0.02 -0.12 0.09 1.001 2921.29 2623.74 0.624
Daily pressure utilized (partner’s view) -0.06 -0.19 0.05 1.002 1774.43 1896.03 0.846
Daily pushing experienced 0.03 -0.02 0.07 1.001 3648.75 3156.17 0.879
Daily pushing utilized (partner’s view) 0.03 -0.01 0.07 1.000 3746.55 2933.40 0.938
Day -0.02 -0.10 0.06 1.000 5367.76 3113.90 0.689
Daily weartime NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 0.07 -0.19 0.36 1.003 635.19 1241.31 0.71
Mean persuasion utilized (partner’s view) 0.04 -0.23 0.32 1.002 651.18 1090.27 0.61
Mean pressure experienced -0.06 -0.38 0.25 1.002 972.96 1406.43 0.655
Mean pressure utilized (partner’s view) 0.00 -0.30 0.30 1.001 1012.72 1369.57 0.506
Mean pushing experienced 0.00 -0.24 0.22 1.007 805.62 1331.43 0.5
Mean pushing utilized (partner’s view) 0.04 -0.20 0.27 1.005 806.08 1396.54 0.643
Mean weartime NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.28 0.21 0.37 1.00 1305.03 2073.99 NA
sd(Daily persuasion experienced) 0.02 0.00 0.06 1.00 1312.22 1681.83 NA
sd(Daily persuasion utilized (partner’s view)) 0.01 0.00 0.04 1.00 2652.14 2196.05 NA
sd(Daily pressure experienced) 0.08 0.00 0.24 1.00 1661.45 2056.83 NA
sd(Daily pressure utilized (partner’s view)) 0.10 0.00 0.28 1.00 1096.63 1593.99 NA
sd(Daily pushing experienced) 0.03 0.00 0.09 1.00 1860.35 1780.90 NA
sd(Daily pushing utilized (partner’s view)) 0.03 0.00 0.09 1.00 1502.56 1673.11 NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA
sigma 0.85 0.82 0.89 1.00 4978.30 2426.63 NA

Reactance

Gaussian

range(df_double$reactance, na.rm = T) 
## [1] 0 5
hist(df_double$reactance, breaks = 6)

hist(log(df_double$reactance + 0.0000001))

formula <- bf(
  reactance ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID), 
  
  hu = ~ persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = ~ cosy(gr = coupleID)
)



prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b"),
  brms::set_prior("normal(0, 2)", class = "Intercept"), # for non-zero PA
  brms::set_prior("normal(20, 10)", class = "Intercept", dpar = 'hu'), # hurdle part
  
  brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0), # Random effects SD for couple
  
  #brms::set_prior("normal(0, 0.3)", class = "ar", lb = -1, ub = 1), # Autocorrelation prior
  #brms::set_prior("normal(0, 0.3)", class = "ma", lb = -1, ub = 1), # moving average prior
  #brms::set_prior("normal(0, 0.3)", class = "cosy", lb = -1, ub = 1), # compound symmetry prior
  
  #brms::set_prior("normal(0.1,2)", class = "sderr", lb = 0),
  brms::set_prior("normal(0.1,2)", class = "sigma", lb = 0) 
)




#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

reactance <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", "reactance_hu")
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(reactance, log_pp_check = T)
## 
## Divergences:
## 5 of 4000 iterations ended with a divergence (0.125%).
## Try increasing 'adapt_delta' to remove the divergences.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 36 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 4000 by 403 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo   -366.1 25.0
## p_loo       108.3  8.7
## looic       732.2 50.0
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.5, 1.6]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     367   91.1%   136     
##    (0.7, 1]   (bad)       33    8.2%   <NA>    
##    (1, Inf)   (very bad)   3    0.7%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(reactance, integer = FALSE)

## Object of Class DHARMa with simulated residuals based on 1000 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help. 
##  
## Scaled residual values: 0.5479685 0.6475247 0.1869277 0.01696357 0.6895541 0.978 0.186306 0.2306857 0.2376249 0.5890393 0.2761681 0.7014718 0.5614676 0.7976995 0.4541524 0.5934053 0.6901254 0.7829586 0.586307 0.7891888 ...
summarize_brms(
  reactance, 
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T,
  include_p_direction = T) %>%
  print_df(rows_to_pack = rows_to_pack_hu)
## Warning: There were 5 divergent transitions after warmup. Increasing
## adapt_delta above 0.8 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Warning in summarize_brms(reactance, model_rows_fixed = model_rows_fixed_hu, :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_direction
Intercept 1.94* 1.49 2.51 1.001 3542.15 3282.10 1
Hurdle Intercept 5.27* 1.89 15.73 1.000 2301.77 3041.25 0.999
Conditional Within-Person Effects
Daily persuasion experienced 0.99 0.88 1.12 1.002 2090.09 2514.53 0.54
Daily persuasion utilized (partner’s view) 1.00 0.89 1.12 1.001 2846.28 2671.13 0.532
Daily pressure experienced 1.00 0.84 1.19 1.003 2732.83 2510.60 0.513
Daily pressure utilized (partner’s view) 1.00 0.65 1.58 1.001 1747.57 1249.07 0.517
Daily pushing experienced 1.03 0.92 1.15 1.001 3747.77 3249.73 0.696
Daily pushing utilized (partner’s view) 1.01 0.82 1.25 1.001 2429.51 2612.71 0.512
Day 1.00 0.66 1.54 1.001 3423.80 3289.05 0.502
Daily weartime NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 0.71 0.41 1.21 1.002 1481.84 2418.50 0.893
Mean persuasion utilized (partner’s view) 1.11 0.66 1.89 1.002 2387.21 2089.54 0.638
Mean pressure experienced 1.55 0.92 2.61 1.001 2173.28 2516.15 0.952
Mean pressure utilized (partner’s view) 1.05 0.63 1.79 1.001 2766.28 2817.06 0.574
Mean pushing experienced 0.90 0.65 1.26 1.000 2072.06 2697.38 0.739
Mean pushing utilized (partner’s view) 1.01 0.70 1.44 1.002 2662.43 2110.30 0.533
Mean weartime NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 1.55* 1.13 2.28 1.000 2495.22 2245.95 0.996
Hu Daily persuasion utilized (partner’s view) 1.03 0.63 1.78 1.002 1848.97 2261.50 0.52
Hu Daily pressure experienced 0.58 0.20 1.89 1.002 2026.14 1798.13 0.869
Hu Daily pressure utilized (partner’s view) 0.96 0.21 6.25 1.001 2407.72 2369.21 0.557
Hu Daily pushing experienced 0.72 0.48 1.05 1.000 2850.30 2636.92 0.955
Hu Daily pushing utilized (partner’s view) 1.18 0.69 2.02 1.001 2827.14 2583.97 0.741
Hu Day 0.75 0.22 2.64 1.002 4534.67 3160.53 0.676
Hu Daily weartime NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 0.63 0.06 4.93 1.006 1166.33 1802.30 0.663
Hu Mean persuasion utilized (partner’s view) 0.71 0.06 6.49 1.005 1264.60 2028.68 0.601
Hu Mean pressure experienced 0.05 0.00 1.21 1.002 1588.17 1929.66 0.963
Hu Mean pressure utilized (partner’s view) 3.45 0.13 159.89 1.002 1673.99 2201.69 0.755
Hu Mean pushing experienced 0.74 0.15 3.69 1.002 1242.70 1850.97 0.66
Hu Mean pushing utilized (partner’s view) 1.22 0.21 7.73 1.001 1421.93 1867.00 0.588
Hu Mean weartime NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.14 0.01 0.36 1.00 1058.25 2219.98 NA
sd(Hurdle Intercept) 1.85 0.99 3.06 1.00 1045.40 1834.72 NA
sd(Daily persuasion experienced) 0.16 0.03 0.31 1.01 881.44 772.35 NA
sd(Daily persuasion utilized (partner’s view)) 0.06 0.00 0.19 1.00 2347.60 2042.87 NA
sd(Daily pressure experienced) 0.10 0.00 0.33 1.00 1903.12 2191.42 NA
sd(Daily pressure utilized (partner’s view)) 0.28 0.01 1.10 1.00 1792.43 1985.13 NA
sd(Daily pushing experienced) 0.09 0.00 0.25 1.00 1125.10 1929.59 NA
sd(Daily pushing utilized (partner’s view)) 0.21 0.03 0.47 1.01 1016.94 815.91 NA
sd(Hu Daily persuasion experienced) 0.47 0.06 0.94 1.00 984.65 983.40 NA
sd(Hu Daily persuasion utilized (partner’s view)) 0.70 0.11 1.46 1.01 1036.50 1088.86 NA
sd(Hu Daily pressure experienced) 1.49 0.18 3.56 1.00 908.71 964.56 NA
sd(Hu Daily pressure utilized (partner’s view)) 1.48 0.07 4.56 1.00 1269.30 1888.72 NA
sd(Hu Daily pushing experienced) 0.43 0.03 1.01 1.00 956.02 1488.17 NA
sd(Hu Daily pushing utilized (partner’s view)) 0.42 0.02 1.17 1.00 1790.38 2251.32 NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA
sigma 0.43 0.35 0.53 1.00 1313.22 1929.56 NA 
hypothesis(reactance, "pressure_self_cw > pushing_self_cw")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (pressure_self_cw... > 0    -0.03      0.12    -0.22     0.17       0.71
##   Post.Prob Star
## 1      0.41     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(reactance, "hu_pressure_self_cw > hu_pushing_self_cw")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (hu_pressure_self... > 0    -0.22      0.63    -1.18      0.8       0.51
##   Post.Prob Star
## 1      0.34     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.

Ordinal

df_double$reactance_ordinal <- factor(df_double$reactance,
                                      levels = 0:5, 
                                      ordered = TRUE)

formula <- bf(
  reactance_ordinal ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = ~ cosy(gr = coupleID)
)



prior1 <- c(
  set_prior("normal(0, 2.5)", class = "b"),
  set_prior("normal(-1, 5)", class = "Intercept"),

  brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0) # Random effects SD for couple
  
  #brms::set_prior("normal(0, 0.3)", class = "ar", lb = -1, ub = 1), # Autocorrelation prior
  #brms::set_prior("normal(0, 0.3)", class = "ma", lb = -1, ub = 1), # moving average prior
  #brms::set_prior("normal(0, 0.3)", class = "cosy", lb = -1, ub = 1), # compound symmetry prior
  
  #brms::set_prior("normal(0.1,2)", class = "sderr", lb = 0)
  #, brms::set_prior("normal(0.1,2)", class = "sigma", lb = 0) 
)




#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

reactance_ordinal <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::cumulative(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", "reactance_ordinal_noar")
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(reactance_ordinal)
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 10 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 4000 by 403 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo   -361.0 24.4
## p_loo        65.8  5.9
## looic       722.0 48.8
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.7, 1.9]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     393   97.5%   195     
##    (0.7, 1]   (bad)       10    2.5%   <NA>    
##    (1, Inf)   (very bad)   0    0.0%   <NA>    
## See help('pareto-k-diagnostic') for details.
summarize_brms(
  reactance_ordinal, 
  model_rows_fixed = model_rows_fixed_ordinal,
  model_rows_random = model_rows_random_ordinal,
  model_rownames_fixed = model_rownames_fixed_ordinal,
  model_rownames_random = model_rownames_random_ordinal,
  exponentiate = T,
  include_p_direction = T) %>%
  print_df(rows_to_pack = rows_to_pack_ordinal)
OR l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_direction
Intercepts
Intercept NA NA NA NA NA NA NA
Intercept[1] 3.93* 1.81 8.79 1.001 4725.90 2920.53 0.999
Intercept[2] 8.00* 3.58 18.08 1.001 4792.84 3156.48 1
Intercept[3] 21.29* 9.04 50.88 1.001 5019.85 3132.24 1
Intercept[4] 90.14* 33.88 262.88 1.001 5516.35 3110.24 1
Intercept[5] 1115.33* 248.00 7355.76 1.000 7530.53 2708.51 1
Within-Person Effects
Daily persuasion experienced 0.67* 0.51 0.85 1.001 4413.62 2680.16 1
Daily persuasion utilized (partner’s view) 1.02 0.72 1.37 1.000 4393.82 2930.89 0.561
Daily pressure experienced 1.57 0.78 2.78 1.000 3564.31 2874.13 0.928
Daily pressure utilized (partner’s view) 1.13 0.43 2.57 1.000 3992.24 2778.44 0.645
Daily pushing experienced 1.31 0.96 1.75 1.000 4946.09 3283.35 0.958
Daily pushing utilized (partner’s view) 0.84 0.57 1.21 1.000 5322.90 2574.58 0.825
Day 1.38 0.49 4.03 1.001 6359.38 3033.27 0.723
Daily weartime NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.24 0.29 5.54 1.001 3064.91 2870.74 0.61
Mean persuasion utilized (partner’s view) 1.26 0.28 6.58 1.002 2847.89 3066.11 0.602
Mean pressure experienced 2.05 0.41 10.86 1.000 4219.38 3194.76 0.819
Mean pressure utilized (partner’s view) 1.01 0.18 5.12 1.000 3986.99 3105.25 0.508
Mean pushing experienced 1.31 0.43 3.83 1.000 3283.62 3364.11 0.7
Mean pushing utilized (partner’s view) 0.77 0.22 2.69 1.000 3362.54 3141.57 0.666
Mean weartime NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 1.15 0.60 1.87 1.00 1755.44 2342.39 NA
sd(Daily persuasion experienced) 0.23 0.01 0.59 1.00 1346.37 1434.84 NA
sd(Daily persuasion utilized (partner’s view)) 0.32 0.02 0.76 1.00 1518.17 2052.41 NA
sd(Daily pressure experienced) 0.76 0.09 1.76 1.00 1530.68 1929.92 NA
sd(Daily pressure utilized (partner’s view)) 0.78 0.03 2.31 1.00 1698.92 2381.24 NA
sd(Daily pushing experienced) 0.36 0.02 0.82 1.00 1061.92 1492.48 NA
sd(Daily pushing utilized (partner’s view)) 0.23 0.01 0.68 1.00 2682.66 2480.40 NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA
sigma NA NA NA NA NA NA NA
disc 1.00 1.00 1.00 NA NA NA NA 
hypothesis(reactance_ordinal, "pressure_self_cw > pushing_self_cw")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (pressure_self_cw... > 0     0.18      0.39    -0.49     0.79       2.23
##   Post.Prob Star
## 1      0.69     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.

Binary

introduce_binary_reactance <- function(data) {
  data$is_reactance <- factor(data$reactance > 0, levels = c(FALSE, TRUE), labels = c(0, 1))
  return(data)
}



df_double <- introduce_binary_reactance(df_double)
if (use_mi) {
  for (i in seq_along(implist)) {
    implist[[i]] <- introduce_binary_reactance(implist[[i]])
  }
}


formula <- bf(
  is_reactance ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = ~ cosy(gr = coupleID)
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b"),
  brms::set_prior("normal(-1, 5)", class = "Intercept", lb=0, ub=5), # range of the outcome scale
  
  brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0) # Random effects SD for couple
  
  #brms::set_prior("normal(0, 0.3)", class = "ar", lb = -1, ub = 1), # Autocorrelation prior
  #brms::set_prior("normal(0, 0.3)", class = "ma", lb = -1, ub = 1), # moving average prior
  #brms::set_prior("normal(0, 0.3)", class = "cosy", lb = -1, ub = 1), # compound symmetry prior
  
  #brms::set_prior("normal(0.1,2)", class = "sderr", lb = 0)
  #,brms::set_prior("normal(0.1,2)", class = "sigma", lb = 0) 
)


#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", "is_reactance_noar")
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(is_reactance)
## 
## Divergences:
## 0 of 4000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 4000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 30 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 4000 by 403 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo   -206.3 12.3
## p_loo        71.0  5.3
## looic       412.7 24.5
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.4, 1.4]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     373   92.6%   102     
##    (0.7, 1]   (bad)       26    6.5%   <NA>    
##    (1, Inf)   (very bad)   4    1.0%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(is_reactance, integer = FALSE)

## Object of Class DHARMa with simulated residuals based on 1000 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help. 
##  
## Scaled residual values: 0.5781259 0.6945909 0.4611728 0.5355948 0.1938734 0.9526177 0.8563463 0.6879258 0.2633464 0.4099011 0.422804 0.07217277 0.3389942 0.08076461 0.4924663 0.3180269 0.02281222 0.3675383 0.0390655 0.03815339 ...
summarize_brms(
  is_reactance, 
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T,
  include_p_direction = T) %>%
  print_df(rows_to_pack = rows_to_pack)
OR l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_direction
Intercept 0.63 0.27 1.47 1.002 3317.92 3032.80 0.86
Within-Person Effects
Daily persuasion experienced 0.62* 0.41 0.86 1.005 2005.26 2215.26 0.999
Daily persuasion utilized (partner’s view) 1.16 0.74 1.95 1.001 2431.93 2416.24 0.732
Daily pressure experienced 1.99 0.69 7.56 1.001 2000.82 1743.82 0.898
Daily pressure utilized (partner’s view) 1.23 0.28 5.87 1.003 2119.85 1169.95 0.627
Daily pushing experienced 1.49* 1.04 2.29 1.003 2079.37 2327.28 0.986
Daily pushing utilized (partner’s view) 0.91 0.53 1.59 1.001 3264.16 2511.94 0.658
Day 1.29 0.41 4.17 1.000 4743.84 3102.51 0.671
Daily weartime NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 3.30 0.49 23.29 1.002 1206.78 1868.46 0.885
Mean persuasion utilized (partner’s view) 2.66 0.38 21.28 1.003 985.82 1731.43 0.823
Mean pressure experienced 14.26 0.82 312.82 1.002 2567.32 2705.68 0.967
Mean pressure utilized (partner’s view) 0.38 0.02 9.35 1.001 2859.03 2720.78 0.73
Mean pushing experienced 1.99 0.41 9.44 1.002 1326.69 2214.51 0.81
Mean pushing utilized (partner’s view) 1.28 0.22 6.79 1.002 1236.69 2699.09 0.616
Mean weartime NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 2.26 1.46 3.32 1.00 1205.19 2002.11 NA
sd(Daily persuasion experienced) 0.46 0.06 0.93 1.00 962.58 1144.68 NA
sd(Daily persuasion utilized (partner’s view)) 0.68 0.10 1.44 1.00 1246.72 1031.43 NA
sd(Daily pressure experienced) 1.39 0.19 3.24 1.00 1008.67 924.37 NA
sd(Daily pressure utilized (partner’s view)) 1.24 0.06 3.52 1.00 1409.06 1769.84 NA
sd(Daily pushing experienced) 0.40 0.02 0.96 1.01 815.40 1446.41 NA
sd(Daily pushing utilized (partner’s view)) 0.41 0.02 1.21 1.00 1321.82 2097.40 NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA
sigma NA NA NA NA NA NA NA 
hypothesis(is_reactance, "pressure_self_cw > pushing_self_cw")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (pressure_self_cw... > 0     0.29      0.65    -0.72     1.43        2.1
##   Post.Prob Star
## 1      0.68     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.

Report All Models

if (report_ordinal) {
  model_rows_random_final <- model_rows_random_ordinal
  model_rows_fixed_final <- model_rows_fixed_ordinal
  model_rownames_fixed_final <- model_rownames_fixed_ordinal
  model_rownames_random_final <- model_rownames_random_ordinal
  rows_to_pack_final <- rows_to_pack_ordinal
} else {
  model_rows_random_final <- model_rows_random_hu
  model_rows_fixed_final <- model_rows_fixed_hu
  model_rownames_fixed_final <- model_rownames_fixed_hu
  model_rownames_random_final <- model_rownames_random_hu
  rows_to_pack_final <- rows_to_pack_hu
}



all_models <- report_side_by_side(
  pa_sub,
  pa_obj_log,
  mood,
  reactance,
  reactance_ordinal,
  is_reactance,
  
  model_rows_random = model_rows_random_final,
  model_rows_fixed = model_rows_fixed_final,
  model_rownames_random = model_rownames_random_final,
  model_rownames_fixed = model_rownames_fixed_final
) 
## [1] "pa_sub"
## Warning: There were 1 divergent transitions after warmup. Increasing
## adapt_delta above 0.8 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Warning in summarize_brms(model, short_version = TRUE, exponentiate =
## exponentiate, : Coefficients were exponentiated. Double check if this was
## intended.
## [1] "pa_obj_log"
## Warning in summarize_brms(model, short_version = TRUE, exponentiate =
## exponentiate, : Coefficients were exponentiated. Double check if this was
## intended.
## [1] "mood"
## [1] "reactance"
## Warning: There were 5 divergent transitions after warmup. Increasing
## adapt_delta above 0.8 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Warning: Coefficients were exponentiated. Double check if this was intended.
## [1] "reactance_ordinal"
## [1] "is_reactance"
# pretty printing
summary_all_models <- all_models %>%
  print_df(rows_to_pack = rows_to_pack_final) %>%
  add_header_above(
    c(" ", "Subjective MVPA Hurdle Lognormal" = 2, 
      "Device-Based MVPA Lognormal" = 2, 
      "Mood Skewnormal" = 2,
      "Reactance Lognormal" = 2, 
      "Reactance Ordinal" = 2,
      "Reactance Dichotome" = 2)
  )

export_xlsx(
  summary_all_models, 
  rows_to_pack = rows_to_pack_final,
  file.path("Output", "AllModels_experimental.xlsx"), 
  merge_option = 'both', 
  simplify_2nd_row = TRUE,
  colwidths = c(38, 7.2, 13.3, 7.2, 13.3, 7.2, 13.3,7.2, 13.3,7.2, 13.3,7.2, 13.3),
  line_above_rows = c(1,2),
  line_below_rows = c(-1)
)
## 
## Attaching package: 'rvest'
## The following object is masked from 'package:readr':
## 
##     guess_encoding
summary_all_models
Subjective MVPA Hurdle Lognormal
Device-Based MVPA Lognormal
Mood Skewnormal
Reactance Lognormal
Reactance Ordinal
Reactance Dichotome
exp(Est.) pa_sub 95% CI pa_sub exp(Est.) pa_obj_log 95% CI pa_obj_log b mood 95% CI mood exp(Est.) reactance 95% CI reactance OR reactance_ordinal 95% CI reactance_ordinal OR is_reactance 95% CI is_reactance
Intercept 55.11* [45.63, 65.99] 121.95* [107.68, 137.86] 4.90* [ 4.80, 5.01] 1.94* [1.49, 2.51] NA NA 0.63 [0.27, 1.47]
Hurdle Intercept 0.34* [ 0.19, 0.62] NA NA NA NA 5.27* [1.89, 15.73] NA NA NA NA
Conditional Within-Person Effects
Daily persuasion experienced 1.03 [ 0.98, 1.09] 1.01 [ 0.98, 1.06] -0.01 [-0.04, 0.01] 0.99 [0.88, 1.12] 0.67* [0.51, 0.85] 0.62* [0.41, 0.86]
Daily persuasion utilized (partner’s view) 1.02 [ 0.97, 1.08] 1.02 [ 0.98, 1.06] 0.00 [-0.03, 0.02] 1.00 [0.89, 1.12] 1.02 [0.72, 1.37] 1.16 [0.74, 1.95]
Daily pressure experienced 0.89 [ 0.77, 1.01] 0.96 [ 0.86, 1.06] -0.02 [-0.12, 0.09] 1.00 [0.84, 1.19] 1.57 [0.78, 2.78] 1.99 [0.69, 7.56]
Daily pressure utilized (partner’s view) 0.92 [ 0.81, 1.03] 0.94 [ 0.86, 1.02] -0.06 [-0.19, 0.05] 1.00 [0.65, 1.58] 1.13 [0.43, 2.57] 1.23 [0.28, 5.87]
Daily pushing experienced 0.98 [ 0.92, 1.04] 1.02 [ 0.96, 1.09] 0.03 [-0.02, 0.07] 1.03 [0.92, 1.15] 1.31 [0.96, 1.75] 1.49* [1.04, 2.29]
Daily pushing utilized (partner’s view) 0.97 [ 0.91, 1.04] 1.01 [ 0.96, 1.07] 0.03 [-0.01, 0.07] 1.01 [0.82, 1.25] 0.84 [0.57, 1.21] 0.91 [0.53, 1.59]
Day 0.92 [ 0.79, 1.07] 0.99 [ 0.88, 1.11] -0.02 [-0.10, 0.06] 1.00 [0.66, 1.54] 1.38 [0.49, 4.03] 1.29 [0.41, 4.17]
Daily weartime NA NA 1.00* [ 1.00, 1.00] NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.08 [ 0.73, 1.58] 1.08 [ 0.78, 1.50] 0.07 [-0.19, 0.36] 0.71 [0.41, 1.21] 1.24 [0.29, 5.54] 3.30 [0.49, 23.29]
Mean persuasion utilized (partner’s view) 1.04 [ 0.70, 1.51] 0.89 [ 0.64, 1.24] 0.04 [-0.23, 0.32] 1.11 [0.66, 1.89] 1.26 [0.28, 6.58] 2.66 [0.38, 21.28]
Mean pressure experienced 0.87 [ 0.39, 1.95] 0.77 [ 0.55, 1.11] -0.06 [-0.38, 0.25] 1.55 [0.92, 2.61] 2.05 [0.41, 10.86] 14.26 [0.82, 312.82]
Mean pressure utilized (partner’s view) 0.71 [ 0.33, 1.58] 0.99 [ 0.71, 1.39] 0.00 [-0.30, 0.30] 1.05 [0.63, 1.79] 1.01 [0.18, 5.12] 0.38 [0.02, 9.35]
Mean pushing experienced 1.21 [ 0.86, 1.70] 1.09 [ 0.84, 1.41] 0.00 [-0.24, 0.22] 0.90 [0.65, 1.26] 1.31 [0.43, 3.83] 1.99 [0.41, 9.44]
Mean pushing utilized (partner’s view) 1.18 [ 0.83, 1.65] 1.17 [ 0.90, 1.52] 0.04 [-0.20, 0.27] 1.01 [0.70, 1.44] 0.77 [0.22, 2.69] 1.28 [0.22, 6.79]
Mean weartime NA NA 1.00 [ 1.00, 1.00] NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 0.65* [ 0.50, 0.81] NA NA NA NA 1.55* [1.13, 2.28] NA NA NA NA
Hu Daily persuasion utilized (partner’s view) 0.75* [ 0.58, 0.92] NA NA NA NA 1.03 [0.63, 1.78] NA NA NA NA
Hu Daily pressure experienced 1.57 [ 0.81, 3.00] NA NA NA NA 0.58 [0.20, 1.89] NA NA NA NA
Hu Daily pressure utilized (partner’s view) 0.69 [ 0.09, 1.98] NA NA NA NA 0.96 [0.21, 6.25] NA NA NA NA
Hu Daily pushing experienced 0.67* [ 0.44, 0.95] NA NA NA NA 0.72 [0.48, 1.05] NA NA NA NA
Hu Daily pushing utilized (partner’s view) 0.58* [ 0.36, 0.84] NA NA NA NA 1.18 [0.69, 2.02] NA NA NA NA
Hu Day 1.28 [ 0.74, 2.20] NA NA NA NA 0.75 [0.22, 2.64] NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 0.68 [ 0.14, 3.48] NA NA NA NA 0.63 [0.06, 4.93] NA NA NA NA
Hu Mean persuasion utilized (partner’s view) 0.69 [ 0.15, 3.45] NA NA NA NA 0.71 [0.06, 6.49] NA NA NA NA
Hu Mean pressure experienced 4.05 [ 0.58, 37.42] NA NA NA NA 0.05 [0.00, 1.21] NA NA NA NA
Hu Mean pressure utilized (partner’s view) 3.31 [ 0.49, 32.67] NA NA NA NA 3.45 [0.13, 159.89] NA NA NA NA
Hu Mean pushing experienced 1.56 [ 0.43, 5.53] NA NA NA NA 0.74 [0.15, 3.69] NA NA NA NA
Hu Mean pushing utilized (partner’s view) 1.13 [ 0.32, 4.21] NA NA NA NA 1.22 [0.21, 7.73] NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.37 [0.28, 0.51] 0.32 [0.25, 0.43] 0.28 [0.21, 0.37] 0.14 [0.01, 0.36] 1.15 [0.60, 1.87] 2.26 [1.46, 3.32]
sd(Hurdle Intercept) 1.47 [1.07, 2.01] NA NA NA NA 1.85 [0.99, 3.06] NA NA NA NA
sd(Daily persuasion experienced) 0.09 [0.04, 0.15] 0.06 [0.02, 0.11] 0.02 [0.00, 0.06] 0.16 [0.03, 0.31] 0.23 [0.01, 0.59] 0.46 [0.06, 0.93]
sd(Daily persuasion utilized (partner’s view)) 0.08 [0.01, 0.14] 0.05 [0.01, 0.11] 0.01 [0.00, 0.04] 0.06 [0.00, 0.19] 0.32 [0.02, 0.76] 0.68 [0.10, 1.44]
sd(Daily pressure experienced) 0.09 [0.00, 0.28] 0.09 [0.00, 0.24] 0.08 [0.00, 0.24] 0.10 [0.00, 0.33] 0.76 [0.09, 1.76] 1.39 [0.19, 3.24]
sd(Daily pressure utilized (partner’s view)) 0.09 [0.00, 0.24] 0.05 [0.00, 0.15] 0.10 [0.00, 0.28] 0.28 [0.01, 1.10] 0.78 [0.03, 2.31] 1.24 [0.06, 3.52]
sd(Daily pushing experienced) 0.05 [0.00, 0.13] 0.09 [0.01, 0.17] 0.03 [0.00, 0.09] 0.09 [0.00, 0.25] 0.36 [0.02, 0.82] 0.40 [0.02, 0.96]
sd(Daily pushing utilized (partner’s view)) 0.07 [0.00, 0.16] 0.04 [0.00, 0.11] 0.03 [0.00, 0.09] 0.21 [0.03, 0.47] 0.23 [0.01, 0.68] 0.41 [0.02, 1.21]
sd(Hu Daily persuasion experienced) 0.26 [0.02, 0.60] NA NA NA NA 0.47 [0.06, 0.94] NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) 0.25 [0.02, 0.56] NA NA NA NA 0.70 [0.11, 1.46] NA NA NA NA
sd(Hu Daily pressure experienced) 0.53 [0.01, 1.76] NA NA NA NA 1.49 [0.18, 3.56] NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) 1.33 [0.08, 3.87] NA NA NA NA 1.48 [0.07, 4.56] NA NA NA NA
sd(Hu Daily pushing experienced) 0.51 [0.06, 1.07] NA NA NA NA 0.43 [0.03, 1.01] NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) 0.55 [0.04, 1.22] NA NA NA NA 0.42 [0.02, 1.17] NA NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA NA NA NA
sigma 0.66 [0.63, 0.69] 0.54 [0.52, 0.56] 0.85 [0.82, 0.89] 0.43 [0.35, 0.53] NA NA NA NA 
report::report_system()

Analyses were conducted using the R Statistical language (version 4.4.1; R Core Team, 2024) on Windows 11 x64 (build 22635)

report::report_packages()
  • beepr (version 2.0; Bååth R, 2024)
  • R.methodsS3 (version 1.8.2; Bengtsson H, 2003)
  • R.oo (version 1.26.0; Bengtsson H, 2003)
  • R.utils (version 2.12.3; Bengtsson H, 2023)
  • brms (version 2.21.0; Bürkner P, 2017)
  • Rcpp (version 1.0.13; Eddelbuettel D et al., 2024)
  • bayesplot (version 1.11.1; Gabry J, Mahr T, 2024)
  • lubridate (version 1.9.3; Grolemund G, Wickham H, 2011)
  • DHARMa (version 0.4.6; Hartig F, 2022)
  • wbCorr (version 0.1.22; Küng P, 2023)
  • tibble (version 3.2.1; Müller K, Wickham H, 2023)
  • R (version 4.4.1; R Core Team, 2024)
  • openxlsx (version 4.2.7.1; Schauberger P, Walker A, 2024)
  • ggplot2 (version 3.5.1; Wickham H, 2016)
  • forcats (version 1.0.0; Wickham H, 2023)
  • stringr (version 1.5.1; Wickham H, 2023)
  • rvest (version 1.0.4; Wickham H, 2024)
  • tidyverse (version 2.0.0; Wickham H et al., 2019)
  • readxl (version 1.4.3; Wickham H, Bryan J, 2023)
  • dplyr (version 1.1.4; Wickham H et al., 2023)
  • purrr (version 1.0.2; Wickham H, Henry L, 2023)
  • readr (version 2.1.5; Wickham H et al., 2024)
  • xml2 (version 1.3.6; Wickham H et al., 2023)
  • tidyr (version 1.3.1; Wickham H et al., 2024)
  • knitr (version 1.48; Xie Y, 2024)
  • kableExtra (version 1.4.0; Zhu H, 2024)
report::cite_packages()